Quantum transfer matrix method for one-dimensional disordered electronic systems 
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We develop a novel quantum transfer matrix method to study thermodynamic properties of one- 
dimensional (ID) disordered electronic systems. It is shown that the partition function can be 
expressed as a product of 2 x 2 local transfer matrices. We demonstrate this method by applying it 
to the ID disordered Anderson model. Thermodynamic quantities of this model are calculated and 
discussed. 
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I. INTRODUCTION 

For real solids, the perfect periodicity is an idealiza- 
tion, while the imperfections are of great importance for 
transport properties. The broken translational symme- 
try makes the system deviate from the extended Bloch 
waves behavior, and in some cases, to localized states. 
As pointed out in Ref. [l[ , we can not adopt the model of 
ordered systems to understand disordered materials. The 
concept of Anderson localization Q and the correlation 
effect Q among electrons in a disordered medium are two 
important ingredients in the understanding of disordered 
systems. 

As a minimal Hamiltonian for independent electrons in 
a disordered potential, the Anderson disordered model 
remains difficult to be understood in finite tempera- 
ture. The disorder itself invalidates conventional analyt- 
ical methods. No proper perturbation parameter can be 
chosen to deal with disordered Hamiltonian although per- 
turbation theory was applied to calculate the conductiv- 
ity in weak disorder limit. Considerable efforts focus on 
Anderson localization and corresponding metal-insulator 
transition. In particular, much insight was gained from 
the scaling analysis [HQ. 

In contrast to higher dimensional cases, ID models are 
often accessible to obtain detailed theoretical (analyti- 
cal and numerical) results. However, in the disordered 
case, it is intractable to make calculation in the ther- 
modynamic limit by analytical methods even in ID. The 
aim of this work is to develop a novel method to re- 
solve this technical problem. Our method avoids direct 
diagonalization of the Hamiltonian and allows the ther- 
modynamic limit to be explored directly and accurately. 
The key point lies in the fact that we can exploit the full 
translational symmetry in the Trotter (imaginary time or 
inverse temperature) direction after trading the evolution 
in real space direction with the Trotter one. 

It should be noted that the transfer matrix introduced 
in the present scheme is not the one usually used in the 
study of disordered systems Q. Our starting point is to 
express the partition function of the system analytically 
in terms of the transfer matrix, rather than to use it to 
trace the eigenvalues or wave functions. 



We will take the disordered Anderson model as an ex- 
ample to demonstrate how the quantum transfer matrix 
method works. The Hamiltonian is defined by 



H = - 5^t i (ctc i+ i+h.c.)+ - aOcJc* 



(1) 



where is the hopping integral between two adjacent 
sites, Ci(cJ) is a fermion annihilation(creation) operator 
at site i, Ui is the diagonal disordered potential, and /j, is 
the chemical potential. U and Ui can take random values 
satisfying some distributions, respectively. 



II. QUANTUM TRANSFER MATRIX METHOD 

As in the quantum transfer matrix renormalization 
group(TMRG) method, we first separates H into 

two parts, H = H\ + H2, with each part a sum of com- 
muting terms: 
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where 



hi.i+i = -ti(cjcj + i + h.c.) + (Ui - [i)m 



(2) 



(3) 



The TMRG uses the second-order approximation of 
the Trotter- Suzuki formula flol [Tl| 



Z = Tr {e- pu ) = Tr (VM)* 1 + 0{e 2 ) 



(4) 



where f3 = 1/fc^T and T is the temperature. (3 is then 
divided into M parts uniformly, e = f3/M and M is the 
Trotter number. 
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where u^i+i are the local evolution operators defined by 



-ehi t i^ 



By inserting 2M identities 
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between the neighboring V\ and V2 operators in (|4|) and 
labeling successively the complete bases with I 6 [1,2M] 
(so called imaginary time's slices), the partition function 
can then be expressed as 
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where V^+l = {n\,n\ +1 \vi^ + i\n\ +1 ,n\\\) represents the 
matrix element of v^i+i. The subscripts i and super- 
scripts I for n and v stand for the coordinates in the 
real and Trotter directions, respectively. If we collect all 
v l 1+1 w hh the same (i, i + 1), the partition function can 
be re-ex pre ssed as the column quantum transfer opera- 
tors 0,1H: 



Z = lim Tr{Ti 2 T 2) 3 • • • TV 1}. 

s— >0 



(8) 



In Eq. ([5]), there exist N site-dependent column transfer 
operators Ti^+i, which are defined by a product of M 
local transfer operators, 
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where the matrix element of the local transfer operator 
t is defined by 



T i,i+1 — 



(ni,l-^ +1 K i+1 |l-^ +1 ,<+i). (10) 
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FIG. 1: (Color online) Graphical representation of the local 
evolution operator v (a) and the local transfer matrix (b) 
along the Trotter and real space directions, respectively. 



The evolution of this matrix is illustrated in Fig. |T]). 
The reason for labeling the basis state by 1 — n\ is to en- 
sures the conservation of the total occupation number of 
two adjacent sites in the Trotter direction. Since h^j+i 
conserves the total occupation number at sites i and i + l, 
subspaces {|1,1)}, {|0,0)} and {|0, 1), |1,0)} are decou- 
pled. In writing Eq. ijTj). a periodic boundary condition 
n\ = nf M+1 is imposed in the Trotter direction. 

Now let uf introduce the following site-dependent vari- 
ables: 
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(11) 



(12) 



Only three of which are independent. In terms of 
these variables, we obtain the following matrix ele- 
ments for Vi i+i in the state number representation 
{|00),|01),|10},|11}}: 

( K 1 \ 

1,1+1 _ , a,i — w l m 
v i,i+l ~°i o Ui a . + w . o 

V oh/ 

According to the procedure mentioned above, the re- 
sulting local transfer matrix in the Trotter direction be- 
tween slices I and I + 1 is also block-diagonal because of 
the fermion number conservation. Consequently, 

( m \ 
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It can be shown that this transfer matrix has the follow- 
ing operator form: 
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where d's are fermion operators defined in the Trotter 
space. T^i is a quadratic function of fermion oper- 
ators. Furthermore, T^i can be exponentiated again 
to a concise quadratic form due to the fermion exclu- 
sion principle. For simplicity, we denote C ; 

Si 



djdi and 



Uib l: 
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then. 
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This exponential quadratic operator form of t ; ;+1 10 
valid only when the four coefficients before (ni — ni + i) 2 , 
(rii — ni + i), djdi+i, and d\ +1 di satisfy the following four 
equations: 



cosh Si 
sinh Si 
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Now we can cast the column transfer matrix T 2 i,2i+i 
defined by Eq. ([9]) into the following form: 



<-i,i+l 



i=i 

+qidl l+1 d 2 i + r l (n 2 i - n 2 i+i) 



(18) 



X^i+i is translational invariant in the Trotter direc- 
tion. Therefore, we can introduce a fourier transforma- 
tion along this direction: 
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The general column transfer matrix can then be rewritten 

as 



T 2 i-i,2i = C£_ x exp 
T 2 i,2i+i = exp 
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Substituting these transfer matrices into Eq. ©, we 
obtain the following expression for the partition function: 



Z = ]JCl*l[Tr[2 + T U) }, 



where 



N/2 
i=l 

t 2 i-i and t 2 i iUJ are 2x2 matrices defined by 
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In obtaining these expressions, we have used the fact that 
T l'l+i/Ci i s an identity matrix in the subspace {|1, 1)}, 
{|0, 0)}. In Eq. (f2"2"]l . all multiplied matrices are expo- 
nential of traceless matrices. Thus the eigenvalue of 



the final matrix after multiplications will have the form 
exp[±A(cj)]. In the thermodynamic limit, the constant 
2 in Eq. (f22|) can be neglected because the eigenvalue 
exp[A(u;)] dominates. 

It should be noted that uj is related to the parity of M. 
For odd M 



M — 1 



.0. 



M — 1 



(24) 
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M 

The transfer matrix T w in Eq. (f2"2"]l is calculated in the 
subspace of {|0, 1), |1, 0)}, so the fermion occupation 
number is 1 for each two unit cells along the Trotter di- 
rection. The state space for one column transfer matrix 
is 



M M , . 

m=-y,...-l,0,...y-l. (25) 



(ft (ft (ft (ft I 



(26) 

and the total occupation number is M. Due to the peri- 
odicity in Trotter direction, it is equivalent to calculate 
in the state space 

44..4 M 4 M+1 \o). (27) 

However, permuting the fermion operator d\ behind d) 2M 
will bring factor (— 1) M_1 because of the transposition 
with [M — 1) occupied fermions. Therefore, there are 
the relation, 

(-1) M " 1 4 = 4 M+1 . (28) 

For odd M, d\ = d\ M+l ; For even M, d\ = —d\ M+l . As 
a result, u> takes values according to the formulas (f24|) 
and ([23)1 respectively. In the following, we will set M to 
even and use Eq. (|2"B"1) . 

As shown in Ref. [13], XL, in Eq. can be formally 
regarded as the transfer matrix for the following ID Ising 
model in a magnetic field: 

Z = 5>xp{^ [Ji(SiS i+1 - 1) + hiSi + fi]}. (29) 

i 

For the odd sites, the parameters are w-independent in 
the corresponding ID Ising model with the relations 



hi-i = --(hi|&2t-i|)) 
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(30) 
(31) 
(32) 



However, the situation is different for even sites and 
the parameters <X, fi are now w-dependent 
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The correspondence for even sites cannot provide a ma- 
terial mapping to ID Ising model because of the complex 
parameters, which arises from the fourier transformation 
along the Trotter direction. Since there exists a relation 
T* = T- u due to the factor e lw , we can multiply in pair 
the transfer matrices with to and — u). In fact, the pair- 
wise multiplication can save us half of the calculation 
because of the complex conjugate relation. When using 
Eq. (|22|) to solve the partition function, what is needed 
is to calculate the product of N matrices of dimension 
2 for a given ui. When Z is obtained, the calculation of 
the free energy is direct, F = —T\nZ (hereafter taking 
k,B = 1), with which all other thermodynamic quantities 
interested such as the specific heat C can be obtained. 

Fig.[2]shows the result of this comparison for the linear 
coefficient of the specific heat C/T. It is obvious that the 
two sets of data coincide completely. The chain length 
TV here and in the rest of the paper is taken up to 2 17 (= 
131072). 




FIG. 3: (Color online) The specific heat C as a function of T 
for a — 0.5, 2, 4, and 0(red solid line). 




o 




FIG. 2: (Color online) Comparison of the specific heat coef- 
ficient C/T obtained from two different methods: the exact 
energy spectrum calculation and the quantum transfer matrix 
method. Here, t = 1 and U = 0. 



III. RESULTS 
A. Gaussian diagonal disorder 

We first consider a Gaussian-like function 

P(x) = -—e^^ (36) 

V27Tf7 

for the disordered distribution of the diagonal potential 
Ui, where o and a are the mean value and the standard 
deviation respectively. In the following discussion, we 



FIG. 4: (Color online) The specific heat coefficient C/T as a 
function of T for a — 0.5, 2, 4, and 0(red solid line). 



take i=l,a = 0, fj, = 0. The controlling parameter 
is the standard deviation a, which denotes the disorder 
degree of the distribution of Ui. In these figures, the red 
solid line means the case without disorder, i.e. Ui = 
(fj = 0), and the other curves are for (Ui) — a = and 

As shown in Figs. [3][5l the difference of the physical 
quantities is slight between the cases of a = 0.5 and 
(j = 0. However, significant differences appear when a in- 
creases. With increasing cr, the peaks of the specific heat 
C move towards higher temperatures. The introduction 
of disordered diagonal energy widens the energy band. 
The density of states(DOS) distributes in a broaden en- 
ergy range. This leads to the shift of the peak position 
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FIG. 5: (Color online) The entropy S as a function of T for 
a = 0.5, 2, 4, and 0(red solid line). 



of C towards higher temperature with increasing a. The 
disorder assists the thermal fluctuations and shifts the 
peak of C/T to lower temperatures. The specific heat 
coefficient C /T at zero temperature is proportional to 
the density of states around the Fermi energy Ep , i.e., 



C 
T 



| T ^o oc p(E F ). 



(37) 



Fig.[4]shows that when the disorder increases, the density 
of state near the fermi surface decreases. The entropy 5* 
decreases when a increases. 

Fig.[5]shows how the chemical potential /i changes with 
T for some given occupation numbers N . N becomes 
large with \i increasing when T keeps invariant, while N 
become small with T increasing when \i keeps intact. To 
keep N invariant, /i is a monotonic increasing function of 
T. The function shows linear shape at high temperature 
regime. 

In Eq. (f2"Tj) . the partition function is expressed as prod- 
uct of different us components, consequently, the free en- 
ergy F can be written as a sum of dependent free energy 
F(u>). For a given temperature T, F(us) decreases with us 
increasing. When the disorder is turned on, F(us) changes 
little except in the vicinity of us ~ tt. This can be seen 
from Fig. [71 which compares the difference of free ener- 
gies between the disordered (er = 2) and ordered cases 
as a function of us. It clearly shows that the difference 
becomes significant only when us approaches tt. In Ref. 
[T3 |. the singularity from some special w's was used to 
discuss the phase transition. 



B. Staggered disorder potential 

We now consider a special model whose diagonal po- 
tential energy is alternating (staggered) with the lattice 



FIG. 6: (Color online) The chemical potential fi as a function 
of T for some fixed occupation numbers N = 0.6,0.7,0.8. 
a = 0, a = 2, and t = 1. 
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FIG. 7: The difference AF(us)/T between disordered and uni- 
form cases as a function of ui/tt. Here, T = 0.01, t = U — 
1, (7 = 2. The inset shows the us dependence of F(us)/T for 
the system without disorder. 



site, i.e. 



k = t, Ui = (-iyu + AUi. 



(38) 



When AUi — 0, the energy spectrum is readily calculated 
via the Fourier transformation, 



p -iki J 



k k 

The result is two sub-bands dispersion relation, 



E± = -2/i ± \/U 2 + 4i 2 cos 2 k . 



(39) 



(40) 



6 




FIG. 8: (Color online) The specific heat C as a function of T 
for the uniform(cr = 0) and disordered(<r = 2, 4) potentials. 
t = l,U = 2,(i = 




FIG. 9: (Color online) Comparison of the specific heat coeffi- 
cient C /T for a continuous distributed random potential with 
that for a discrete random potential, t = 1,U — 2, /j, = and 
T = 0.01. 



The band gap is 2U. 

Let us choose a uniform rectangular distribution 
P(x) = — for the increment of the random diagonal en- 
ergy AE/j. Here, a is the width of the rectangular distri- 
bution satisfying 



if 



P(x) 



— - < X < - 

9 — x — 9 1 



0. 



(41) 



otherwise. 



In the presence of disorder, the band gap is expected to 
decrease with the increase of the disorder degree a. For 
small a, there is a finite excitation gap and the specific 
heat drops exponentially in low temperature, as shown 
in Fig. M ■ 

We also consider the case when the random poten- 
tial take only two discrete values: —a/2, a/2 with equal 
probability. We call it discrete distribution of the incre- 
ment, and the above rectangular distribution is denoted 
by continuous one. 

Fig. [5] compares the specific heat coefficient C/T at 
T = 0.01 for the above two kinds of distribution of ran- 
dom potentials. In the case of discrete distribution, C/T 
exhibits a sharper peak. When the disordered level a 
increases to approximately 4, the band gap disappears. 
With further increasing a, C drops to zero because the 
band gap open again for the discrete random potential. 
This can be understood as follows. The two sub-bands 
close to each other when disorder is introduced. The up- 
per sub-band shifts downwards by a/2 and the lower sub- 
band shifts upwards by a/2. When the top of the original 
upper sub-band touches the zero energy, i.e., a/2 = v6, 
the two sub-bands begin to separate again. Therefore, 
C/T decreases to zero at a ~ 2y6 in Fig. [9l On the con- 
trary, for the continuous random increment case, C /T 
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FIG. 10: (Color online) Temperature dependence of the 
specific heat in three different random potentials(<j = 
3.0,4.8,7.0). t = l,U = 2,n = 0. 



remains a finite value even for large a. This is because 
the splitting of the upper and lower bands only blur- 
ringly expands the width of these two bands, and once 
they touch each other, they never separate again. 

As a further investigation for C in discrete distribu- 
tion, we choose three typical disordered degrees a = 
3.0,4.8,7.0, corresponding to three regions in Fig. [51 to 
show how C varies with the temperature at low T. As 
shown in Fig. [TQl when a = 3.0 and 7.0, C drops ex- 
ponentially with temperature in low temperatures. The 
double-peak structure of C can be understood from the 
overlaps of energy bands. The discrete increments split 
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the original single energy band into the two bands which 
shift upwards and downwards by a/2. The resulting four 
bands from upper and lower sub-bands meet pairwise. 
The thermal fluctuations bring double peaks shown in C 
curves. 

IV. CONCLUSION 

For the ID disordered system, the quantum transfer 
matrix method we have developed is applicable to all 
kinds of disorder distribution types and strengths. The 
non-diagonal disordered problem can be handled since 
the partition function Z can be expressed as the product 
of site-dependent local transfer matrices. Compared to 
the diagonal disordered cases, we only need to modify 
the local transfer matrix elements correspondingly. 

We have studied the thermodynamic properties of the 
ID disordered Anderson model. We discussed two kinds 
of diagonal (potential) disordered models with or with- 
out staggered potentials. The free energy F can be writ- 
ten as a sum of different lo components from the Fourier 
transformation in Trotter space. Comparing with the 
system without disorder, the most significant difference 



in F(us) shows only in the region very close to to — n. The 
disorder change the distribution of DOS, leading to the 
difference in the thermodynamic quantities in compari- 
son with the disorder-free system. All the results shown 
in the figures are for systems with the number of sites 
greater than fO 5 . This kind of calculation is far beyond 
the capacity of exact diagonalization. 

The transfer matrix method has a broad range of appli- 
cability and can be used to discuss any non-interacting 
fcrmion models. Recently this method has been used 
to calculate the thermodynamic quantities of the Hof- 
stadter model which describes the behavior of tightly- 
bound Bloch electrons in the mag netic field 0,[l5|. In 
Landau gauge, the Hofstadter Hamiltonian can be decou- 
pled into a sum of one dimensional Hamiltonian, which 
lies in the application range of the quantum transfer ma- 
trix method. 
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